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I. INTRODUCTION 



The discovery of high temperature superconductivity in the copper based transition metal 
oxides has revived the interest in superconductivity of strongly correlated electron systems. 
Mechanisms of high temperature superconductivity which are not phonon-mediated have 
been sought for a long time, and their existence has remained a controversial subject. This 
old problem cannot be settled by perturbative methods 

The large-dimensionality limit of the strong correlation problem has received much re- 



cent attention §,[g,|[,[g,§,0|[[|ig[y] especially in the context of the Hubbard 
model. This is a nontrivial limit of the strong correlation problem which can be solved 
non-perturbatively. 

In this paper we define a model with two degrees of freedom per unit cell which also 
has a well defined limit when the coordination number gets large. On the atomic scale the 
interactions in this model are purely repulsive. We show that the normal phase is stable 
against phase separation and displays interesting metal-insulator transitions as a function 
of doping, and of the interaction strength. 

For the first time we give in this a paper a derivation of the self consistent equations 
for the local Green's functions in the presence of off-diagonal superconducting long range 
order that become exact in the limit of infinite dimensions. We then solve these equations 



by extending the Hirsch-Fye Quantum Monte Carlo algorithm |12 



We present evidence that in a well defined region of parameters of the two-band model 
the normal state is unstable against a superconducting state, characterized by a strongly 
frequency dependent superconducting order parameter. We observe that superconductivity 
in this model is promoted in regions of parameters where two different atomic configurations 
are close in energy, that is in the mixed valence regime. 



II. THE MODEL 



We consider the following hamiltonian: 



[dfj^'Pa) represent two atomic orbitals on different sublattices (A, -B) of a bipartite lattice. 
The ('copper') orbital d„ is strongly correlated, while the ('oxygen') orbital Po- is uncor- 
related. Each site has identical connectivity so that the model describes a 'CwO'-type 
system. A two-dimensional three-band model has been proposed by Emery [l^ and by 



Varma et al [18| as a minimal model of the copper-oxide Cu02 planes. The model we 
consider is similar in spirit, while having a non-trivial large d limit (as opposed to other 
possible generalizations of multiband models which have somewhat degenerate large d limits 
0, [|19|)- Multiband models have been investigated intensively, but no clear numerical evi- 



dence for superconductivity in two dimensions has so far been obtained ||20|, In infinite 
dimensions several circumstances facilitate the numerical work: the model reduces to an 
effective single-site model, in which the thermodynamic limit is taken exactly, without the 
need for finite-size extrapolations, and the Monte Carlo algorithm is free from minus-sign 



problems, allowing us to reach lower temperatures than in previous works ||23[| . 
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In the absence of correlations {Ud = 0), diagonalization of Ti yields two bands (bonding 
and antibonding) : = {ep + e^i w (cp — e^)^ + 4e^}/2, where ek is the Fourier transform of 



tij. These bands are separated by a gap Aq = ep — e^. (Having in mind the copper-oxides in 
the hole representation, we choose throughout this paper < Cp). For a total hole density 
n = J2a{.< ^da > + < npa- >) Smaller than 2, the Ud = ground state is a partially filled 
copper band and an empty oxygen band. For 2 < n < 4, the Ud = ground-state is a full 
d-band and a partially filled p-band. The T = chemical potential is discontinuous as a 
function of density at n = 2, with fi{n = 2~) = ed, fi{n = 2"*") = e^. Hence, for Ud = 0, this 
model describes a metal for all densities except n = 2, where one has a band insulator. The 
copper and oxygen density of states have simple expressions for Ud = 0: 



Ndie) = f-^ Ni^ie - ep)(e - e,)) iVp(e) = X-^ iV(^(e - ep)(e - e,)) (2) 
V e — * \^ e — tp ' 

where N{e) = J2k6{e- ek). 

The limit of infinite connectivity z ^ oo requires a scaling of the hybridization tij as: 
Uj = tpd/y/z, SO that the density of states A^(e) has a proper limit 0. In practice, one may 
consider the d-dimensional hypercubic lattice {z = 2d), for which: A^(e) = l/-\/27rtp^ 
as (i oo or the Bethe lattice with connectivity z for which A^(e) = — {e/tpdY /2T[tpd as 
z oo. We shall establish the general equations for an arbitrary N{e) but, for the sake of 
simplicity, will perform all numerical simulations for the Bethe lattice. 

It is by now well-established that, in the limit ^ oo, the problem reduces to a single-site 
'impurity' model supplemented by a self-consistency condition. Standard methods that 
will not be reviewed here allow us to derive the corresponding equations for the one-particle 
Green's functions. In this section, we assume no long-range order of any kind (magnetic or 
superconducting) and display equations valid in the paramagnetic normal state. 

Since only the 'copper' sites are correlated, all the local quantities can be derived from 
an impurity model for those sites only. The effective action of the impurity model is given 
by: 

S = Ud f dTnd^{T)ndi{T)- f dr dr' Y.d^{r)D,\r - r') d^r) (3) 
Jo Jo Jo ^ 

We denote by D{r — r') = — < Td{T)d^{T') >s the interacting Green's function calculated 
with this action and by the impurity self-energy (seen as a functional of Dq), T^diioJn) = 
DQ^{iujn) — D~^{iujn)- The self-consistency equation for Dq then reads, {un = (2n+ l)7r//?): 

r N(e) 

D{iuJn) = Cp(^Wn) / de ;/ (4) 

where we have set set (p^iujn) = ioJn + — Cp, Cd(^'^n) = i^n + — ^d — ^d{i^n)- 

Once Eqs.(H) ,(^ are solved for Dq , the impurity self energy evaluated at the self 
consistent value of Dq gives the d electron lattice self energy In the ((iko-jPko-) basis, it reads 
in matrix form: 



1 



CpCd - 




(5) 
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In the z oo hmit, self-energies become purely site-diagonal, so that depends only on 
frequency and the Sp^ component is absent. The absence of the diagonal component Sp 
in eq.(^ comes from the simplifying assumption of an uncorrelated p-orbital. From eq.(^), 
one sees that the self-consistency equation simply means that the impurity model Green's 
function must coincide with the on-site d-orbital Green's function: D{iun) = J2]iD{k,iujn)- 
Also, the p-orbital on-site Green's function is simply given by: 

Pitur,) = O(^^n) / de f (6) 

It is straightforward to extend these equations to a model involving a Hubbard repulsion 
on oxygen orbitals also. One then has to solve simultaneously two separate impurity models, 
one for each orbital, and both and Sp enter the self-consistency equations (§J^). It is 
also possible to include a direct oxygen-oxygen hopping tpp (with a.l/z scaling as z ^ oo). 
For the sake of simplicity, these additional terms will not be considered in this paper. 

Throughout the rest of this paper, we shall work with the z = oo Bethe lattice with the 
semi-circular d.o.s. given above. This is motivated by the very simple form taken by the 
self-consistency equations in this case: 

= iujn + fx-ed- tl^P{iuJn) P'^ = + - ep - tl^D^iu^) (7) 

A direct way to derive these equations for the Bethe lattice is to integrate out fermionic 
degrees of freedom on all lattice sites except a single one. Because of the z = oo limit, 
only one-particle processes (involving the full Green's functions D, P) are generated in this 
partial summation, yielding the effective action eq. (0), with Dq explicitly given in terms of 
D by eq.(0). 

In order to solve the coupled problem defined by eqs.(^, (^, we use a numerical method 
that has been described in detail elsewhere [0] |[10| ||rT|. It is based on an iterative 
procedure: for a given imaginary-time Do{t), the impurity model (^) is solved (given a 
discretization of the interval [0, jS] in L = (3 /At time slices) using the Hirsch-Fye quantum 
Monte-Carlo algorithm. This produces a discretized D{t) which, after Fourier transforma- 
tion, is inserted into eq.(|^ in order to produce a 'new' function Dq. This process is iterated 
until convergence is reached. 



III. METAL-INSULATOR TRANSITIONS 
A. The metal-insulator transition at n = 1 

For a density n = 1 this model displays a transition from a metallic to an insulating 



state following a scenario outlined by Zaanen, Sawatzky and Allen . 

For very large Ud/Ao and tpd « Aq = ep — e^, the n = 1 ground state has all copper 
sites nearly singly occupied and all oxygen sites nearly empty. Hybridization between the 
two orbitals costs an energy ~ Aq, while the kinetic energy gain is only ^ tp^/Ag. Hence, 
for tpd << Ao, we have a charge-transfer insulator at n = 1, with a jump in the chemical 
potential fi{n = 1^) — ^{n = 1^) ^ 0. The numerical data for nd,np and the total density 
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n as a function of the chemical potential fl = fi — ed confirms these expectations. In fig. |l|, 
we show data at /5 = 30, Ud = 8, and Aq = 4. 

In the opposite limit (tpd » Aq) the delocahzation energy wins, and we have a strongly 
correlated metal with a strong hybridization of 'copper' and 'oxygen' orbitals. A metal to 
charge-transfer insulator transition separates these two regimes, at (Ao/tpd)c = 0(1) for 
large Ud- 

In the opposite regime of weak correlations {Ud << Aq), the metal-insulator transition 
at n = 1 has the character of a Mott transition within the copper band. Assuming that 
tpd « Aq, most of the d-orbital density of states is, for Ud = 0, concentrated around e^, with 
a small bandwith of order tp^/Ao (see eq.^ and the hybridization with the oxygen orbital 
is weak. As Ud is increased, the lower copper band is gradually split by the interactions, 
and a Mott transition occurs when Ud becomes comparable to the bandwith tp^/Ao, i.e. 
(Ao/tpd)c — tpd/Ud for small Ud- In fig. ^ we show data at /? = 30, Ud = 1.5, and Aq = 4 
displaying this behavior. 

Based on these arguments, the phase boundary separating the metallic and insulating 
regimes at n = 1 is expected to have the schematic shape |]13[ described in fig. ^. 

When f/ ~ Aq, there is a crossover regime which to our knowledge has not been investi- 
gated previously. Data in this region is shown in fig. ^. The detailed investigation of these 
phase transitions is left for future work. 



B. The metal to band insulator transition at n = 2 

Interactions can also induce a transition from insulating to metallic behavior at a density 
n = 2. At Ud = the system is a band insulator with a band gap Aq = ep — e^. Turning on 
the interactions has the effect of reducing the band gap. Whenever the gap is less than the 
temperature, we expect metallic behavior. 

This change from insulating to metallic behavior at n = 2 is clear in the nvsfi curves in 
figs |l] and |. In fig.0 the case Ud < Aq is realized, for which the atomic ground-state |T4 



at n = 2 has nearly all copper sites doubly occupied and all oxygen sites empty. We see gaps 
at n = 1 and n = 2 with fi{n = 1"*") — fj,{n = 1~) ~ Ud, /i(n = 2+) — yu(n = 2~) ~ Aq — Ud- 
In fig ID the case Ud > Aq is realized and the corresponding atomic ground-state at n = 2 
has all copper and oxygen sites singly occupied. Here, gaps exist at = 1 and n = 3, 
with /i(n = 1+) — /i(n = 1") ~ Aq, /i(n = 3+) — fi{n = 3") c::^ Ud — Aq. 

Hence an insulator to metal transition can be thought to occur at n = 2 at a critical 
value of the ratio Ud/ Aq ~ 1, which can be interpreted as the upper Hubbard band hitting 
the p level. In fig. (^ such an intermediate situation is realized. Below, we find that 
very interesting physics arises in the mixed valence regime {Ud — Aq), where the strongly 
correlated copper level lies close to the uncorrelated oxygen level. 

Additional insights on this transition can be gained by looking at the self-consistency 
equation (^ for the special point /i = ep, which reads for low-frequency {iun — >■ cj-t-iO"*", u — > 
0): 

D(0 + zO+)-i + Sd(0 + zO+) = ep-ed + D{<d + zO+)-^ (8) 

The low-frequency behaviour of D and Ud is severely constrained by this equation, which 
allows two different solutions: i) D(zO+) = and ii) Srf(zO"'") + = ep. Case i) is obviously 
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realized for small Ud, while case ii) can be interpreted as the effective d level E(iO^) + 
becoming degenerate with the p level. Although we cannot obtain the numerical solution at 
zero temperature we find numerical evidence that one goes from i) to ii) as Ud is increased. 
This is readily seen in the data of fig.^- § for four different values of the interaction strength. 
From these curves, it is apparent that T,d{iO~^) remains close to — for all values of Ud 
larger than a critical value. These figures display in fact the results of two subsequent 
iterations of the self consistency loop, which we display in order to highlight the remarkable 
accuracy of the Monte Carlo algorithm. 

The qualitative behavior of n vs /i can also be understood from a Hubbard Ill-type 
solution of the coupled equations (0), which takes into account the physics of high-energy 
processes in a qualitatively correct manner [P], [0. The appropriate generalization of the 



Hubbard HI approach to the two-band model inserts the following expression for D{iuJn) 

D{tun) ^ + J } (9) 

Z Uq — Ud 

into the self-consistency equation (|^), to obtain the following closed equation for D{iuJn)'- 

- Uj/A) + tl,{xpUj/2 - 1xpx\ + tl^Xd) 
+x^ - Ujxl/A + tl^XpXd + tprfx) D-xXp = Q ^ ' 

where Xd = iuJn - ed - Ud/2, Xp = iun - ep, x = XpXd - t^^. 

This simple approximation gives a rather satisfactory account of the various metal- 
insulator transitions described above: zero-temperature spectral densities are displayed in 
fig ^ for various values of the parameters. It is unable however to account for the low-energy 
quasiparticle excitations in the metallic regimes, since it just ignores the Kondo effect in the 
impurity model eq. d^). 



IV. SUPERCONDUCTING PROPERTIES 

A. Dynamical pairing in infinite dimensions 

In this section we extend the d = oo formalism to incorporate the possibility of supercon- 
ducting long-range order. The equations that we will establish here allow the investigation 
of the model within the broken symmetry phase ( cf where the Hubbard model was 



studied within magnetically ordered states). Let us define 'anomalous' Green's functions: 

Fdir - t') = -T < d,,i(r)d,T(r') > F,(r - r') ^ -T < p,i(r)p,,t(r') > (11) 

describing the formation of on-site pairs. Singlet pairing corresponds to F even: F{t) = 
F{—t) = —F{l3 — t), while Sz = triplet pairing corresponds to F odd: F{r) = — F(— r) = 
F{i3 — r). Allowing for a non-trivial time-dependence of F is crucial. The underlying 
physical idea is that on-site equal-time pairing is likely to be strongly suppressed in the 
presence of a strong on-site repulsion but that pairing involving a 'time-lag' between the 
paired holes may occur. This idea dates back to Berezinskii's proposal |]21]] for triplet pairing 



in ^He, a generalization of which has been recently considered for cuprate superconductors 
by Balatsky and Abrahams in the singlet case. 
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In the presence of a non-zero F it is convenient to work with Nambu spinors = {d^ , rfj 
(similarly \&p) and with the matrix formulation of one-particle Green's functions: 



D(. - r') . -T < Mr^ir') > = ^^Z^, ^,ZJ. (12) 



With these notations, the kinetic term of the hamiltonian Ti reads: —tij'^^^a^'^ p,j where (73 
denotes the Pauli matrix. Following the usual method, we integrate out fermionic variables 
on all sites except on a single copper site. The 'impurity' action obtained in this way now 
reads: 

Ssup = Ud 1^ dTnd^{T)ndi{r) - dr T rfr'^+(r)Do-^(r - r')^d(r') (13) 
Jo Jo Jo 

where Do is given in terms of D and P by the self-consistency equations: 

Bo'^iiuJn) = i^n + (/i - ed)0-3 - (T3P(^C^n)0^3 
V~^{iuJn) = iuJn + (/U - ep)0-3 - t^^(T3D(it^„)(T3 

We can account for an externally applied dynamic pairing field A(i(zco'„) on all copper sites 
in the original lattice problem by adding a forcing term 



Ad(iu;„ 
Arf(it^„) 



(15) 



to the r. h. s. of eq. (|T^). 

The impurity action eq. ( pTSf ) describes an Anderson impurity in a superconducting 
medium. Since this problem, even with static pairing turns out to be highly non trivial, 
we can expect that the self consistent solution of eq. (|T^ will allow very intricate densities 
of states. 



B. Algorithm in the presence of pairing [16] 



In the present paragraph we show how the Hirsch-Fye algorithm |]T2| can be generalized 
in the presence of off-diagonal terms (in the spin indices) in the Green's function eq. (|12D 
which, after the usual Trotter breakup, becomes a 2L x 2L matrix. 

The quartic term in the action is decoupled, as usual, by the introduction of Ising vari- 
ables. 

exp[— Arf/n^n-'' — [ri' + n^)/2 = Tro- exp[A(T(/)[n^ — n^\\ (16) 



According to the prescription given by Soper |[T5|, all terms in the action have to be 
normal-ordered (here, in ^f"*", \l/ with \l/'^(r) = (\I'j['"(r), '^ti'^)) = (cf, Cj)), before performing 
the discretization. Consequently, eq. (|16D is written as 

Tr<,exp[A(T[^+^i + ^+^2] exp[-A(T(/)] (17) 

The Green's functions for different spin configurations are related by the following Dyson 
equation: 
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D,, = + (D, - l){exp{W' -W)- 1)D,, 



(18) 



where Da-(r, r') is the full 2L x 2L matrix Green's function of \l/ for a given spin configuration 
(j(/), and Wk = Wk+L = \cr{k). Eq. (p!8| ) is derived in the same way as the one in ref. |p!2 



and differs from that expression only in a sign for the down-spin sector. The Green's function 
of the model is given by 



D.(r, r')Det{-D~') exp(-A Ei ^(/)) 
E./^et(D-i)exp(-AEiO^(0) 



(19) 



The statistical weight of a configuration is thus given by the product of the usual fermion 
determinant, and of a scalar factor, which arises from the commutator in eq. ([TB|). 

As in the Hirsch-Fye algorithm ||12| , we can use the Dyson equation to calculate the ratio 
R of statistical weights for flipping a spin cr(A;) cr'{k) = —a{k), 



R = R^Riexp{-2Xa[k]) 



(20) 



R^ = l + [1-Bik,k)]f (21) 

R^ = 1 + - -D{k + L,k + L) - D(fc + L, fc)D(fc, k + L)f/{1 + (1 - D(A;, k))f)]f (22) 

with / = [exp(— 2Acr(A;)) — 1]. Once a flip of spin a{k) is accepted, all the elements of the 
2L X 2L Green's function are updated in a 2-step procedure, corresponding to the flip of Wk 
and Wk+L- 

D'(Z, n = D(/, /') + (D(/, k) - 6(1, k))f/{l + (1 - D(A;, fc))/)D(fc, /') (23) 

D"(/, /') = D'(/, /') + D'(/, k + L)- 6{l, k + L))f/{1 + (1 - B'{k + L,k + L))f)B'{k + L, I') 

(24) 

We have extensively compared the algorithm with exact diagonalization results for an 
impurity interacting with a few conduction electrons in the presence of explicit pairing terms 
d'^dtfj, and conclude that the algorithm is correct, and practical. In fig. [l^ we show data 
of one such realization of the superconducting Anderson model for At = 1, .5, .25, .125, and 
(exact diagonalization). Notice the quadratic convergence. 



C. Numerical Results 

To investigate the transition into the superconducting state at finite temperature we 
define a pairing susceptibility Xp{i^n,i^m) 

Fd{iuJn) = '^Xp{^^n,il^m)Ad{iuJm) Ad{iUm) (25) 
m 

This susceptibility is a well defined symmetric matrix and is finite in the normal state. 
The divergence of its largest eigenvalue signals the transition into the superconducting state. 
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In our approach it is natural to express Xp iii terms of two operators, and a. XxiiuJn, i^m) 
describes the pairing response of the impurity model eq. ( [I3l) 



Fd{iuJn) =^\K{ii^n,i^^rn)[GQd ^]l2(«t^m) [God ^]l2^0 (26) 
m 

while the diagonal operator a describes the response of the self consistent medium to the 
impurity site and can be obtained from eq.(0) as a(n, m) = 6n,mtpdP*{'i^n)P{i^n)- 
Clearly, Xk, Oi and Xp ^-^e related by the following equation 

Xp=[I~\KaY'\K (27) 

When the largest eigenvalue of A = Xxd approaches one, the susceptibility diverges. It is 
clearly more accurate to study the behavior of the largest eigenvalue of A = AfcO? and we 
do that by the power method using the full algorithm. To do this, we run the algorithm 
described in the previous section a large number of times, starting with a small value of 
the off-diagonal term -F(r, r') of the Green's function. In the normal state, -F(r, r') will 
converge to zero and two subsequent runs will yield two functions functions F(r, r')and 
F'(r, r') obeying 

F'(r, t') = J dT"\{T, t")F{t", t') ~ Xma.F{T, t') (28) 



Xmax is the largest eigenvalue of the matrix A (c/eq. ( [26| ) and F'(r — r') converges to 
the corresponding right eigenvector, which is also the eigenvector corresponding to the most 
unstable eigenvalue of the susceptibility matrix. 

In fig (|ll|), we show Xmax vs n at a relatively high temperature = 15 in the singlet 
sector. I.e., after explicitly anti-symmetrizing the functions F{t) and Fq{t) with respect 
to (3/2. In order to assess the importance of finite Ar effects, we have plotted Xmax for 
L = 16, 32, 64. Let us stress that this extrapolation is the only error in our calculation, since 
the statistical uncertaintly and the uncertainty stemming from the self-consistency condition 
are extremely small. All of our data is consistent with quadratic convergence in Ar^ for 
small values of Ar which we have clearly established for the Green's function itself. 

The susceptibility seems to have a local maximum close to the density n = 1, and is 
again large at around n = 2. We notice that the tendency towards superconductivity is 
largest near half filling and near ~ 2, i.e. in the mixed valence regime. 

We have been able to produce perfectly converged superconducting solutions of the mean- 
field equations in the singlet case for temperatures within our reach. The solution for F(r), 
at temperatures just below Tc(Ar) is very similar to the largest eigenvector of the suscepti- 
bility matrix above Tc. The finite value of Ar produces large finite size effects in the value 
of Tc(Ar). At this point we have established that Pc{L = 64) = 37, Pc{L = 128) = 59 and 
f3c{L = 192) = 75. A reliable finite size scaling of this quantity to determine the value of the 
critical temperature requires larger simulations and will be described elsewhere. Here we 
use the fact that the largest eigenvector of the susceptibility matrix converges very well to 
discuss the nature of the possible superconductiviy in this system. Of particular importance 
is the different and strong frequency dependence on the copper and the oxygen sites. 

Examples of such numerical solutions are given in fig |1^ for 128 time slices, at the value 
which seems to be most favorable for the appearance of superconductivity Aq = 4. The 
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Fourier transform of the order parameter on the copper site near Tc are displayed in fig 



13 . Notice (inset of fig 12) that the pairing amphtude is much larger on oxygen sites than 
on copper sites, as expected, and that the equal-time pairing -F(O) is nearly zero on copper 
sites. 



V. CONCLUSION 



In this paper we introduced a controlled approach to the study of the occurrence of su- 
perconductivity in strongly correlated systems. We extendended the mean field approach to 
the strong correlation problem of ref. 0] to derive a self consistent picture of a correlated su- 
perconductor which becomes exact in infinite dimensions. As in the corresponding treatment 
of the normal state ( |^) an impurity model plays the role of the mean field hamiltonian. 
The impurity model corresponding to a correlated superconductor is an Anderson impurity 
model embedded in a superconducting medium. The quantities that correspond most closely 
to the Weiss field is now a set of two functions that describe the diagonal and the off diagonal 
response of the effective medium. A numerical solution of the self consistent equations for 
the Green's functions revealed a strong frequency dependence of the superconducting order 
parameter. 

Using this method we investigated the most promising regions for high temperature 
superconductivity in the present version of the extended Hubbard model in large dimensions. 
The largest tendencies towards superconductivity were found when the normal phase is is 
in a mixed valence regime. 

There are several natural extensions of this work. One can study the effect of adding 
extra interactions in the hamiltonian on the superconducting transition temperature of the 
system. While the critical temperatures we found are relatively high, a better understanding 
of the behaviour of the solution of our equations may lead to higher transition temperatures. 
Finally, as emphasized throughout the text, larger simulations are called for. They should 
allow the finite size scaling analysis which can fully elucidate the existence of the super- 
conducting state proposed here. A completely open problem is the study of the physical 
properties of this superconducting state. We stress that because of the strong correlations a 
non perturbative treatment like the one outlined in our paper is necessary for tackling this 
problem. 
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Figure Captions 

1. Hp (full line) rid (dotted line) and n (dashed line) V5 /i. — at Aq = 4 [/^ = 8, and 

2. Up (full line) (dotted line) and n (dashed line) vs ii — aX = A = 1-5, and 
/3 = 30. 

3. The Zaanen Sawatzky Allen phase diagram . The special point occurs when Ud ~ Aq 

4. Hp (full line) na (dotted line) and n (dashed line) — at Aq = 4 C/^ = 4.5, and 

= 30. 

5. Im{G{iuJn)) (full line) and Re{Tid{i'^n)) (dotted hne) for Ud — 1-5, — = 4, Aq = 4, 
and/3 = 45 (L = 128). 

6. Im{G{iuin)) (full line) and Re{T^d{i^n)) (dotted line) for Ud = 2.5, — = 4, Aq = 4, 
and/3 = 45 (L = 128). 

7. Im{G{iu}n)) (full line) and i?e(Sd(ia;„)) (dotted line) for Ud = 4.5, — = 4, Aq = 4, 
and/3 = 45 (L - 128). 

8. Im{G{iuJn)) (full line) and i?e(S(i(ia;„)) (dotted line) for Ud — 6.5, — = 4, Aq = 4, 
and = 45 (L = 128). 

9. Hubbardlll approximation for the one-particle density of states Pd(a;) vs uj for Aq = 4, 
and Ud = .25,3,4.5,8 

10. Green's function G{t = 0)and F{t = 0) in one realization of the superconducting 
Anderson model for At = 1, .5, .25, .125, and (exact diagonalization). We have 
observed quadratic convergence of the Green's functions for any value of r and a wide 
range of parameters. 

11. X vs n at P = 15 and L = 16 (upper curve), 32, and 64 (lower curve) at the special 
point Ao = 4:,Ud = 4.5, - = 4, 

12. Green's functions Gdij) and Fd{T) at the special point Aq — A, ji — ed — ^, at Ud — 4.5. 
/? = 64 and L = 128. The inset shows Fd{T) and Fp^r). 

13. The real part of Fd{iuJn)vs ujn on the special point at L = 128 and inverse temperatures 
13 = 64(full hne), (3 = 70(dotted line), and /? = 80(dashed line), 

14. Imaginary part of Gd{ioJri) on the special point at L = 128 and inverse temperatures 
13 = 64(full hne), (3 = 70(dotted line), and (3 = 80(dashed line). 
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